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We present a differential rate equation model of chiral polymerization based on a simple copoly- 
merization scheme in which the enantiomers are added to, or removed from, the homochiral or 
heterochiral chains (reversible stepwise isodesmic growth or dissociation). The model is set up for 
closed systems and takes into account the corresponding thermodynamic constraints implied by the 
reversible monomer attachments, while obeying a constant mass constraint. In its simplest form, 
the model depends on a single variable rate constant, the maximum chain length N, and the ini- 
tial concentrations. We have fit the model to the experimental data from the Rehovot group on 
lattice-controlled chiral amplification of oligopeptides. We find in all the chemical systems employed 
except for one, that the model fits the measured relative abundances of the oligopetides with higher 
degrees of correlation than from a purely random polymerization process. 
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I. INTRODUCTION 



In the transition from prebiotic racemic chemistry to chiral biology one scenario suggests that homochiral peptides 
must have appeared before the appearance of the primeval enzymes [H 0] ■ While several stochastic synthetic routes 
for mirror symmetry breaking that convert racemates into nonracemates have been described [3, 4] , the generation of 
long bio-like polymers [lj made up of repeating units of the same handedness requires elaboration of new synthetic 
routes. Polymerization reactions of racemic mixtures of monomers in solution are typically expected to yield polymers 
composed of random sequences of the left- and right-handed repeat units following a binomial or Bernoulli distribution. 
Thus the probability for obtaining oligomers with homochiral sequence becomes negligible with increasinglength [l| . 

Recent investigations have proposed that TV-carboxyanhydride (NCA) @ and thioester derivatives 0,13] of amino 
acids might have operated as relevant precursors for the formation of the early peptides [8£. Results on the polymer- 
ization of NCA monomers in organic solvents, [ij-fhl] in water [l5l-[i~7| and in the solid state [HI, Oil have been published. 
Luisi and coworkers [2(j| - |23j have reported the polymerization of racemic a-amino acids in solution which yields small 
amounts of oligopeptides of homochiral sequence whose abundances with respect to the heterochiral chains exhibit a 
slight departure from the binomial distribution. 

This problem of the random distribution can be overcome by catalyzed polymerization of amphiphilic amino acids, in 
racemic and nonracemic forms, which self-assemble into two-dimensional ordered crystallites at the air-water interface 
[24l . Based on a process involving self-assembly followed by lattice controlled polymerization, Lahav and coworkers 
recently proposed a general scenario for the generation of homochiral oligopeptides of a single handedness from non- 
racemic mixtures of activated alpha amino acids [241 [25j . Initial non-racemic mixtures undergo a phase separation 
by self-assembly into a 2D racemic crystalline phase and a separate enantiomorphous 2D phase of the enantiomer 
in excess. Each of these crystalline phases has markedly different chemical properties, thus yielding products that 
differ in the composition of the oligomers. So, polymerization within the enantiomorphous crystalline phase yields 
homochiral oligopeptides of one handedness whereas the reaction controlled by the racemic crystallites yields racemic 
mixtures and heterochiral products. The combination of the two routes leads to an overall chiral amplification process. 

In this paper, we are interested in the lattice-controlled polymerization reactions proposed by those authors. It 
is important to clarify at the outset what specific aspect of the overall experimental mechanism we want model 
here and the way we aim to do so. The proposed experimental scheme starts from an initial excess, say S > R 
of monomers which undergoes an initial self-assembly process into two types of two-dimensional crystallites at the 
air/water interface. Once formed, each one of these two crystal phases participates in the control of a subsequent 
type of polymerization. Thus, the racemic crystallites polymerize racemic mixtures of oligomers and the heterochiral 
products, whereas the other pure enantiomorphous crystallite controls the polymerization of the isotactic chains, these 
are formed from the monomer in excess (S, in this example). However, the details of the polymerizations depend in a 
complicated way upon the specific packing arrangements of the crystal monomers and the possible reaction pathways 
taken within each crystallite phase. The authors of the experiments state that the connection between the monomer 
packing arrangements in the crystallites and the resultant composition of the various diastereoisomeric products is 
"not straighforward" (25[. We therefore opt for a simple model for interpreting their data. With this objective in 
mind, we present a copolymerization model for the interpretation of the experimental data. The model may be 
termed effective in the following sense: it presupposes or takes as given the prior formation of the self-assembled 2D 
crystallites at the air-water interface and is concerned exclusively with the subsequent polymerization reactions. Thus 
the complicated microscopic details referring to the monomer packing arrangements and reaction pathways within 
the crystallite self-assemblies are treated implicitly with our rate constants. Our copolymerization reaction rates can 
satisfactorily account for the different chemical properties of the two crystalline phases (racemic 2D crystallites and 
pure enantiomorphous 2D crystallites) that lead to the formation of racemic mixtures, heterochiral products and 
isotactic oligopeptides. We contrast the fits from our model with those assuming a purely random process that obeys 
a binomial distribution. The final justification for considering such an effective model rests on its ability to yield 
good fits to the data. The goodness of the fits obtained below demonstrates that the experimental data can be fit 
convincingly as if the simple scheme depicted pictorially in Fig. Q]were the sole mechanism leading to the observed 
relative abundances. This then gives additional meaning to to term "effective" , and in the operational sense. 



Our starting point is a simple model for the copolymerization of two chemically distinct monomers displaying a 
wide variety of product sequence compositions. The model we introduce and study here is an appropriately modified 
and extended version of the one considered a few years ago by Wattis and Coveney [26| . 

The main important differences compared to prior and related models are that we (1) consider polymerization in 
closed systems [271]-, so that no matter flow is permitted with an external environment- and (2) we allow for reversible 
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THE COPOLYMERIZATION MODEL 




3 



* Homopolymerization 




• Heterodimerization 




• Heterotrimerization 




• Heteropolymerization 




FIG. 1. The copolymerization model. The (R)-chiral (red) and (S)-chiral (blue) monomers reversibly associate into the growing 
homochiral (top) or heterochiral (bottom) copolymer chains. Because the system is closed, both the heterodimer (second line) 
and hetero-trimer (third and fourth lines) reactions must be treated separately to avoid double counting and thus ensure that 
the total system mass is conserved in a closed system (see text for an explanation). 



monomer association steps. We also correctly include the formation (and dissociation) of the heterodimer [27[. It 
turns out this must be treated on a separate basis in order to avoid double counting, which if left unchecked, would 
lead to a violation in the constant mass constraint. Once the heterodimer is treated correctly, this implies that the 
hetero-trimer must also be treated separately. Beyond this, the remainder of the hetero-oligomers can be treated in 
a uniform way. 

First, we introduce the notation to be used. Polymers are classified by three quantities: the number of A monomers 
of which it is composed (subscript r), the number of B monomers which it contains (subscript s) and the final or 
terminal monomer in the chain, denoted by a superscript. In this scheme, the monomers are denoted by A = C^ 
and B — Cq X ; pure homopolymers are denoted by C^ and Cq s ; all copolymer chains C^ s or C^ s with r, s > 1 are 
heteropolymers. Note also that chains of the form Cq b and C^ are forbidden. The corresponding time-dependent 
concentrations are denoted by lower case variables: e.g., c^ s (t) and c^ s (t). The model is then defined by the following 
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reactions, in which equilibrium is maintained between the finite monomer pool and the ensemble of copolymers: 

k 

C r . s + A ^ - C r+1 , (1) 

k* 

-, \ kab i , 

C r,s + B ' C r,s+H ( 2 ) 



■6 



C£+B^^C£ +1 . (4) 



This model can accommodate any two chemically distinct monomers. For the purpose of this paper, we consider the 
case when A = R and B = S are two enantiomers. 

The overall basic scheme must be broken down into several special subcases, especially important so as to avoid 
undesired double counting of the heterodimer and heterotrimer reactions, see Fig. [1] Once we treat these special 
cases, we then pass to the corresponding set of rate equations for the concentrations. 

The formation of chirally pure polymer chains denoted by c^ and Cq u , for 1 < n < N — 1 is described by the 
homo-polymerization reactions: 



riA , n A aa t (~<A n B ■ s~tB bb „ n B 
°n,0 ' °1,0 °ra+1.0 °0,n "+" °0,1 ^ °0,n- 

V ^bb 



(5) 



iV is the maximum chain length permitted. In our recently reported work [27jj , we considered that once a monomer 
has been added to a homopolymer of the opposite chirality (that is, "the wrong" monomer), the polymer is inhibited 
and further growth is halted. This polymer could not directly react anymore and could only lose its wrong terminal 
monomer through the inverse reaction. In the present model, we assume such a chain can continue to grow by adding 
monomers of both configurations. So, for 2 < n < N — 1, the hetero-polymerization or inhibition reactions are as 
follows: 



f~u\ i r ,ts - c a r< a _i_ r< A - /~<a 

u n,0 °0,1 ^ " u n,l u 0,n "+" °1,0 ^ " 

Kb Ka 



(6) 



For both homo- and hetero-polymerization reactions, represented by Eq. (5JH1 the upper limits specified for n ensure 
that the maximum length for all oligomers produced (or consumed) by these reaction sets, both the homo- and 
heterochiral ones, is never greater than N. In the remainder of this paper we will consider here the natural and chiral 
symmetric reaction rate assignments k aa = kbb, k ao — k oa and likewise for the inverse rates, k* a = k bb and k* b = k ba , 
reducing the number of independent rate constants to four. 

Even if we have the information about the composition, we can only know the chirality of the last monomer attached 
to the chain, we have no information regarding the specific sequence. This implies that the following two reactions 
are indistinguishable: 

a D k a b D d 4 kba a 

CA i r~iB ^ f^B r^B i r~iA w 
i.o "r °o.i °i.i °o.i + °i.o 

k* k* 
^ab K ba 

Thus for all practical purposes, = C-f^ and this suggests using the following notation: C\ \ = Cf x = C-f^ and 

to define a unique direct constant rate: kh = kah + kha , and an inverse one A£ = ha . Note that if k a b — kb a , then 
kh = k a f, = kta- Due to these characteristics, we will treat the heterodimer in a different way compared with the other 
hetero-polymers. The reaction of the heterodimer formation is therefore: 

kh 

Ci,o + Qi.i - - ^1,1 • ( 7 ) 

K 
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As before, the reactives and products in Eq. ([7]) are the same, so the differences in the free energy between initial 
and final states should be the same in all the reactions in these equations, implying the following thermodynamic 
constraint on the reaction rates: 



^ab kba /q\ 

pr - tt- ( 8 ) 

K ab Ha 

If the heterodimer formation were not to be treated in the separate way as we have done, and were to be included, 
e.g. in Eq.© by merely changing the lower limits for n (2 < n < N — 1) by 1 < n < N — 1, we would be making the 
mistake of double counting it. The same occurs for the heteropolymers formed from the addition of a monomer to a 
heterodimer. The two reactions of each pair of the following equations are also indistinguishable: 

A k aa A A kba A 
1.1 +W0 T °2 1 l^l.l tOjqt L/ 2 i 

K aa Ha 



g kbb g £> ^afc g 

Ci,i + C o l ^ =^ C l2 Ci,i + c a i ^ c± 



Hb ft at> 



Again, it is convenient to define the following direct reaction rates for these steps, kha 



k aa +k ba 7, _ k bb +k a 
9 i ^hb — 9 



and inverse k^ a = fc "°+ fc t,a , k^ b — kbh ^ kab . Note that if k aa = k ob and k ao = ki, a , then k^ a — k) l0 , and if k* a = k% b and 
Kib = Ka' then ^ha = ^hb- The reactions to consider are then: 

A ^ha a q kfoh ^ 

Cl,l + Cl ' — ^2 1 Cl.l +^01 ' —^12- (9) 

Ha K hb 

As we have already remarked, in our model, as in the original one for open systems |2q |. the polymeric chains that 
have taken up the "wrong" chirality monomer can continue to grow. Thus, we we allow for the further growth of 
these chains by adding monomers of either chirality. This kind of polymerization reaction for 2 < n < N ~ 2 is given 
by: 

a a k aa . B B b 

C\, n + C 10 ■ > — C 2 C n l + C 01 v - C n 2 (10) 

K aa Hb 



C A + C B 
°l,n + °0,1 



r B 



C B 



°1.0 



(11) 



And for2<r<iV-2,l<s<iV-l-r: 



^r,s T b l,0 : 



ft A 



■a 



0,1 



Hb 



r< B 

■ u r,s+l 



(12) 



C 



A 

r.s 



k a b 



■ °r,s+l 



a 



B 
r.s 



/~iA 
°1,0 ' 



kba 



(13) 



Note that in the elementary reaction steps, in the rate constants, and in the corresponding differential rate equations 
(see below), the left-right symmetry of the model is manifest, that is, possesses a discrete Z 2 symmetry. This symmetry 
can be broken spontaneously by the dynamical solutions of the differential rate equations, thus this model is apt for 
studying spontaneous mirror symmetry breaking. 

By lifting the Z 2 degeneracy in the reaction rates, e.g., allowing for k aa ^ kbb and thus leading to more independent 
rate constants for describing the reaction set, we could study the influence of explicit chiral bias in the model. As this 
is not the aim of this work, we will not consider it here. 
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We next write down the differential rate equations corresponding to this reaction network. We employ the rate- 
equation theory as in chemical kinetics. We begin with the rate equations for the chiral monomers: 

, A N-l N-2 N-2N-l-r N-l N-2 N-2N-1-S 

ac l,0 



— fcgqc^o ( 2ci + y^ o + y^ c^„ + y^ y^ fc6 a c^ f y^ + y^ 4 + y^ y^ 

n=2 n=2 r=2 s=l n=2 n=2 s=2 r=l 

AT jV-2 N-lN-r 

— khCi QCy^ — kha.CipCi.1 + k* a (2C2 + E C^ o + ^ + ^ ^ C^ Jj 



n— 3 n— 2 r— 3 s— 1 

AT-1 JV-2 W-lW-r 

*6a (Ei+E C t+EE 4>) + fc ^ C l.l + k ha4l ( 14 ) 

n=2 n=2 r=3 s=l 



JV-1 AT-2 W-2Af-l-s AT-1 N-2 N-2N-l-r 

"'< u -k bbC ^ Oic^ + e c ln + E c «,i + E E <&) - ( E <o + E + E E < 



n=2 n=2 s=2 r=l n=2 n=2 r=2 

N N-2 N-\N-s 



dt 

— khCi oCn i - khbCQ^ci^i + kl b (2cQ 2 + E + E c ™> 2 + E E c ^ s ) 

n— 3 n— 2 s— 3 r— 1 

JV-1 AT-2 JV-lJV-s 

+ fc: b ( E <i + E <2 + E E + fc ^i.i + 2 (is) 



n=2 n— 2 s— 3 r— 1 



The equations describing the concentration of the homopolymers, for 2 < n < JV — 1: 

— = hbCoA^C^n-l ~ C 0%) — k ba CQ n cf Q + ^C^„ +1 — C^„^ + fc£ a C^„ (17) 

It is necessary to treat the kinetic equations of the maximum length homopolymers N individually. Since these do 
not elongate further, they can not directly react, and can not be the product of an inverse reaction involving a longer 
chain: 

^ — k a aCifiC N _ 10 — k* a c N0 (18) 

dc 0N 

^ ~ kbbC . 1 c Q N _ 1 — kl b c Q N (19) 

The differential equations describing the concentration of each type of heteropolymer (included the heterodimer) , for 
2 < n < N - 2 : 

^ = kh c tfl c o,i ~ k ha ci.ic^ Q — khbCi^CQ^ — k* h cx.\ + k^c^i + k^cf^ (20) 



dr A 

ac l,n 

dt 



— —k aa c 10 c ln — k ab c 01 c ln + k ba c 0n c 10 + k* aa c 2 , n + k* b c ln+1 k^ a c ln (21) 

= — ^M> C 0^1 C ra,l — kbaCi fiCn,! + ^ab c n,0 C M k bb c n,2 + ^ha c n+l,l — ^afc c n,l (22) 



As before, it is useful to treat individually the maximum length polymers iV: 



dc A 



- /c 6o c 0,JV-l c l,0 K ba c l,N-l 



dt 

. A B i * B 

K a6CjV-l,0 C 0,l — K ab C N _ 11 



dr B 

ac N-l,l 



dt 

As was mentioned when describing the reaction network, each kind of trimer c A A and c B 2 must have its 
differential equation in terms of kh a , khb- 



dc A A 
dt 

dc B 2 
dt 

For 2 < n < N - 3 : 



'k a aC A fi C 2,l kabCQiC^ + k} la C\^C A + k aa C A x + k ab C 22 k ha C 21 

-kbbC B iC B 2 — k ba c A Q c B 2 + khbCi t \c B i + k bb c B 3 + k ba c 22 — k hb c B 2 



da 



dt 
dt 



r~ ~ k aa c A (c A n c 2,nj k a bC B iC A n + kb a c B n c A + k aa (c A n c 2,nj + k ab c 2n+1 k ba c 2n 
kbbC a i(c nl — c n2 ^j — k ba c 10 c n2 + k ab c nl c 01 + k bb (c n3 — c„ ;2 ) + k ba c n+12 — k ab c n2 



Once again, the equations corresponding to the maximum length homopolymers N are 

dr A 



dt 



dt 

For 3 < r < iV - 2 and Ks<N-l-r 



— k aa CifiC 1N _2 + k ba c 1N _ 2 c 10 k^ a c 2 N _ 2 — K a c 2 N _ 2 

dr B 

ac N-2,2 _ , B B _i_ U U( A „B _ u* B _ , * B 

ft -bfc c 0,l C N-2,1 ' ft ab c 7V-2,l c O,l Hb c N-2,2 K ab C N-2,2 



dc A 



j£ — kaaCiQ (c r _i s c r,s) k ab C 01 C rs +k ba C r _ ls C 10 + k aa ^C r+ls C r ,s^j + ^ab c r,s+l k ba C rs 



dt 



For 3 < s < N - 2 and f < r < N - 1 - s: 



dc B 



dt 

For 3 < n < N - 1 : 



— kbbCQA (^ c r,s-l C rys) ^ba c lflC B s + k ab C A s _ 1 C B 1 + k bb (c B s+1 C r,s^ + ^6o C H-l,s ^ab c r,: 



dr A 

UL -n,N-n _ , a A , . R „A ; * A 



,,. ^ao C l,0 C n-l,AT-n + ^&a c n-l,N-n c l : ^aa c n,N-n ^ba c n,N-n 

i 

J-n,n _ i B B . i A B i* B ,* B 

— — K &6 c 0,l c iV-n,7i-l T K abC 7v _ n rl _ 1 C ^ Kbb c N-n,n K ab c N~n,n 
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TABLE I. Number of differential equations as a function of the maximum polymer length N 





Number of eqs 




Number of eqs 


c l,0 


1 


c 0,l 


1 


Cl,l 


1 






Cn,o, (2 < n < N) 


^ ' N ATI 

£„=2 = JV - 1 


co, n , (2 < n < N) 


E„ =2 = N - 1 


Ci, n , (2 < n < N - 1) 


E„= 2 = W - 2 


c„,i, (2 < n < N - 1) 


E„=2 =AT-2 


c£,i 


1 


Cl,2 


1 


C2,n, (2<n< iV-2) 


^ — ^iV — 2 7i r o 

E„= 2 = N - 3 


Cn,2, (2 < n < N - 2) 


^ — \N — 2 TV r o 

E„= 2 =AT-3 






c B 




(3 < r < N - 2) 


Ef= 3 2 Ef^ = i(iV 2 - 7iV + 12) 


(3 < s < N - 2) 


E^- 1 1 - s Ef= 3 2 = K^ 2 -7iv + i2) 


(l<s<jV-l-r-) 




(l<r<JV-l-s) 




<N-n, (3<n<JV-l) 


Et-^iv-3 


c?,s, (3 < n < N - 1) 





As remarked earlier, the complete reaction scheme must satisfy mass conservation in a closed system, implying that 
the mass variation rate must be strictly zero: 

N N-l 

= 2c M + 3(c^ + cf 2 ) + ]T n« + c* n ) + + !)(^n + <£,i) 

JV-2 JV-1JV-1 

+ £ (n + 2) (c^„ + c£ 2 ) + ^^(r + S ) « s + c* s ) , (36) 

n— 2 r— 3 s— 1 

where the overdot stands for the time-derivative. The compliance with this constraint is an important and crucial 
check on the consistency of the numerical integration of the full set of differential equations Eqs. (I14M35[) . which 
we monitor and confirm in all the simulations presented below. Analytically, this relation is satisfied by the rate 
equations. 

As we see, there is one differential equation for each type of monomer and one for the hctcrodimcr. The homopolymer 
set requires E/n=2 — 1) equations and the heteropolymer set a total of 2(N — 2) equations. The total number 
of kinetic differential equations describing the whole system is N(N + 1), and is broken down into the separate 
contributions as displayed in Table HI Then, the total number of equations for describing the system as a function of 
maximum chain length N is: 

#eqs = 6 + 2(A^ - 1) + 2{N - 2) + 2{N - 3) + (N 2 -IN + 12) + 2{N - 3) = N(N + 1), 

(37) 

as pointed out in Ref [27j . From the computational point of view, the number of equations grows quadratically with 
the maximum chain length N. 

III. NUMERICAL RESULTS 

We are interested in applying our copolymerization model to fit the experimental data measured by the Rehovot 
group, so our primary goal is to reproduce as closely as possible the details reported concerning the experiments on 
chiral amplification of oligopeptides. For this purpose, the first step is to determine the initial monomer concentrations 
to be employed in the simulations. The actual experiments were carried out for 0.5mM solutions of monomers, thus 
we have employed for each case: (a) R : S = 1 : 1 which corresponds to an initial enantiomeric excess eeo = 0%, so 
c^ (0) = 0.25toM and cg A (Q) = 0.25toM; (b) R : S = 4 : 6 corresponding to ee = 20%, so c^ (0) = 0.2mM and 
c£l(0) = 0.3mM; (c) R : S = 3 : 7 which corresponds to ee = 40%, so cf o (0) = 0.15mAf and 0^(0) = 0.35mM. 
The remainder of the initial concentrations (the dimers and on up) are taken to be zero. Next, we systematically 
search for the reaction rates leading to the best fit to the given data. 

Different chemical model systems were used in the experiments: namely 7-stearyl-glutamic thioethyl ester (Cis — 
TE — Glu), iV e -stearoyl- lysine thioethyl ester (Ci$ — TE — Lys), 7-stearyl-glutamic acid N-carboxyanhydride (Cig — 
Glu — NCA) and 7-stearyl-glutamic thioacid (Cis — thio — Glu), varying both their initial compositions and for 
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various choices of catalyst. The composition of the oligopeptides formed was analyzed by matrix- assisted laser des- 
orption/ionization time-of-flight mass spectroscopy (MALDI-TOF) with enantio-labeled samples. The experimental 
relative abundances of the oligopeptides was inferred from the ion intensity. It are these relative abundances that we 
aim to interpret vis-a-vis our copolymerization model. 

Since only the experiments with racemic mixtures of the starting compounds required a catalyst, it is reasonable to 
expect that the racemic and the chiral enriched cases will follow different dynamics for a given model system. That is, 
the presence or absence of a specific catalyst affects the rate constants, for a given chemical system. Firstly, we will 
find the reaction rates for the racemic case, and afterwards, those for the enriched chiral case, allowing us to compare 
both. The a-priori nine free parameters we must set to run the numerical integrations are comprised by the four direct 
and the four inverse rate constants k aa , fc^, k a b, kb a , an d k* a , k bb , k* b , k ba , plus the maximum polymer chain length, 
N. We set all the inverse reaction rates to a unique value, fc* a = k bb = fc* b = k ba = 10 _10 (s _1 ), implying an almost 
irreversible scheme, and we determine the remainder of the parameters from fitting the copolymerization model to the 
relative abundance data. This required numerical integration of the set of differential equations Eqs. ([14II35P which 
we performed using the Mathcmatica program package. For each independent run we verified the compliance of the 
numerical results with the constraint in Eq. ([36[l . an imperative for any closed system. 

Results from fitting the model to the data indicate that the maximum chain length TV does not play a significant 
role, the Pearson product-moment correlation coefficient, r, remains the same for N = 12, 14, 16, 18, 20, so we will set 
N = 12 for all compounds and cases treated below. Since the number of independent equations scales as ./V 2 , this 
represents an important reduction on computer time and the memory used. We note that one is free to scale out the 
dependence of one pair of reaction constants from the rate equations by a suitable redefinition of the time variable. 
Thus, without loss of generality, we set the cross inhibition rates equal to unity k a b = kb a = l(s^ 1 mol~ 1 ) and then 
search for the reaction rates k aa = kbb leading to the best fits. 

A. Racemic mixtures 

In one set of experiments, the authors reported MALDI-TOF analysis of the oligopeptides formed at the air-water 
interface from racemic mixtures R : S = 1 : 1 of the monomers for the various model systems and catalysts. We first 
fit the copolymerization model to this data. 

The best correlation data for the racemic Cis-TE-Glu system, with the I2/KI catalyst are found for k aa = kbb = 
1.7(s -1 moZ -1 ). In this case, the best fit obtains for the time scale t = 10 11 (s). Exactly by the same process, the 
best correlation data for the racemic C\s — TE — Lys are found for k aa = kbb — 2.3(s~ 1 mol~ 1 ) and for k aa = 
kbb = 1.3(s -1 moZ ) when adding 1%/KI and AgNOz as catalyst, respectively. For the simulations here, we took 
the times t = 10 10 (s) and t = 10 11 (s) in the racemic cases with I2/KI and AgNO^ respectively. Finally, we fit our 
copolymerization model to the Cig — thio — Glu experimental relative abundances. The authors of the experiments 
affirmed that this compound undergoes a truly random polymerization, so fits from our model are expected to be 
slightly less satisfactory than those for the binomial distribution function. Setting the inverse reaction rates and the 
cross inhibition as indicated above, then the best correlation coefficients are found for k aa = kbb = 0A(s mol ). 
The instant or time-scale leading to these numerical values is t = 10 10 (s). 

The corresponding (experimental and numerical) relative abundances for the four compounds cited above corre- 
sponding to these values are shown in Fig. [21 The histograms show the relative abundance of each experimentally 
obtained oligopeptide compared to the best fit from our copolymerization model. We emphasize that we fit the model 
to the complete family of stereoisomer subgroups (global fit). The resulting data correlations are shown in Fig. [3] and 
Table [Til the latter gives a detailed comparison of the best fits between individual subfamilies and the overall global 
fit. 

In the case of the Cig — Glu — NCA with catalyst Ni{CH^C02)2, the best fit is obtained for k aa = kbb = 
0.2(s~ 1 mol~ 1 ). Results for the corresponding relative abundances are shown in Fig. [Hand the correlation from fitting 
is displayed in the bottom frame of Fig[3] and Table IIIII Not all subfamily data sets are reported in the experimental 
paper [241 ] ; here we use the fitted model to the partial data set to predict or fill in this missing subfamily data. 
Numerical results for the racemic case have been found for t = 10 10 (s). 



B. Chirally enriched mixtures 

In a second set of experiments, the authors reported MALDI-TOF analysis of the oligopeptides formed at the air- 
water interface from non-racemic mixtures of the monomers for the same model systems. No catalysts were employed 
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RS C,o-TE-Glu 1:1 



Copolymerization Model 




FIG. 2. Relative abundance versus number of repeat units (r, s) of the oligopeptides obtained from fitting the model (white) 
to the experimental data (black) from racemic mixtures R : S = 1 : 1 of monomers. The four chemical models are indicated by 
the insets. 



TABLE II. Comparative fits between the copolymerization model and the binomial distribution to the experimental relative 
abundances: racemic mixtures R : S = 1 : 1 of monomers of the four model systems as indicated in the leftmost column. Only 
in the case of Cis — thio — Glu does the binomial distribution give a better global fit than the copolymerization model: this 
latter system provides an experimental reference system for random polymerization [25l ]. 



r 


Copolymerization model 


Bin. 


Fits for each subgroup n 


Global 
fit 


Global 
fit 


di 


tri 


tetra 


penta 


hexa 


Cis -TE- Glu 


0.92 


0.96 


0.80 


0.84 




0.93 




Cig -TE- Lys(h/KI) 


0.96 


-0.82 


-0.11 


-0.73 


0.45 


0.85 


0.32 


Cis — TE — Lys(Ag) 


0.98 


1 


0.03 


0.88 


0.76 


0.84 


0.8 


Cis — thio — Glu 


1 


1 


1 


0.98 


0.97 


0.95 


0.98 



there. We next consider fits of our model to these data sets. 

The best correlations factors for both chirally enriched mixture cases (20% and 40% excesses) in the case of the 
Cis — TE — Glu system are found for the same rates, that is for k aa = kbb = 2(s~ 1 mol~ 1 ). The results for these 
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Model 



FIG. 3. Data correlations r from fitting the model to the data in the case of racemic mixtures of all the compounds employed. 
The chemical systems are indicated by the insets. The solid line represents the linear correlation between experimental data 
and numerical calculations. 



values are shown in Table IIVI In Fig. [5] we display the relative abundances of the homochiral oligopeptides and in 
Table \V\ both the calculated and experimental enantiomeric excesses for the 4:6 and 3:7 (R:S) mixtures. In Fig. [6]we 
show the data correlation. Numerical results for the non-racemic case have been found for the time scale t — 10 11 (s). 

For the chiral mixtures of C\g — TE — Lys we found the best fits for the dynamics corresponding to k aa = fc^ = 
2.5 (s^mol^ 1 ). The results for these values are shown in Table IVT1 The relative abundances results for these values 
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RS C 18 -Glu-NCA 1:1 
Copolymerization Model 



A I 



I 



FIG. 4. The Cig — GZm — NCA system with catalyst Ni(C HsCOz)^. relative abundance of the oligopeptides obtained from 
fitting the model (white) to the experimental data (black) from racemic mixtures of monomers. Compare to Fig. 4A of 
reference [2^] 



TABLE III. Comparative fits between the copolymerization model and binomial distribution to the experimental relative 
abundances for the racemic compositions (R:S=1:1) of the Cis — Glu — NCA system. 



r 


Copolymerization model 


Bin. 


Fits for each subgroup n 


Global 
fit 


Global 
fit 


di 


tri 


tetra 


penta 


hexa 


hepta 


octa 


nona 


deca 


endeca 


dodeca 


Cig - Glu - NCA 


1 




1 




0.98 




0.98 




0.97 


1 


0.95 


0.96 


0.75 



TABLE IV. Comparative fits between the copolymerization model and the binomial distribution to the experimental relative 
abundances measured for non-racemic mixtures of Ci$ — TE — Glu. 



r 


Copolymerization model 


Binomial 


Fits for each subgroup n 


Global 
fit 


Global 
fit 


di 


tri 


tetra 


penta 


hexa 


(R:S) 4:6 


0.86 


0.89 


0.93 


0.99 




0.94 


0.75 


(R:S) 3:7 


0.95 


0.94 


0.96 


0.99 


0.99 


0.95 





TABLE V. Enantiomeric excesses ee: numerical results from the copolymerization model (experimental data) for the relative 
abundances of the homochiral oligopeptides for the Cig — TE — Glu system. 



ee(%) 


di 


tri 


tetra 


penta 


hexa 


(R:S) 4:6 


18 (26) 


24 (39) 


30 (46) 


35 (59) 




(R:S) 3:7 


37 (48) 


48 (71) 


57 (82) 


66 (92) 


73 (>99.8) 



are shown in Fig. [7]and the enantiomeric excesses obtained for 4:6 and 3:7 (R:S) mixtures are presented in Table IVIll 
In Fig. |S] the data correlation is shown. For the simulations here, we took the instants t = 10 10 (s) and t = 10 11 (s) in 
the racemic cases with I2/KI and AgNOs respectively, and t = 10 10 (s) for the chirally enriched mixtures. 

In the case of nonracemic C\$ — thio — Glu, the best correlation coefficients are found for the same values of the 
reaction rates that we found in the racemic case, namely for k aa = = 0A(s~~ 1 mol~ 1 ). Results for the chiral cases 
are shown in Table ["VTlTl As to be expected and as shown there, the correlation factors for the global fit to the binomial 
distribution function are slighter better than those for any simulation we could perform with the copolymerization 
model, so we reconfirm what was claimed by the authors of the experimental work: namely that the Cis — thio — Glu 
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FIG. 5. Relative abundance versus number of repeat units (r, s) of the oligopeptides obtained from fitting the model (white) 
to the experimental data (black) from non-racemic mixtures of monomers for the Cis — TE — Glu system. 



Data Correlation 



Q. 
X 

w 0.1 
0.6 
0.4 



1 1 1 I 1 1 

T ■ B 


■ C18-TE-Glu4:6 

y-« 

y= 1.015X 

r = 0.94 


^— i 1 1 1 1 1 1 

■ • • 

■_ i i i 


1 J-l 

■ C18-TE-GIU 3:7 

y = x 

y= 1.068 

r= 0.95 



0.0 0.1 0.2 0.3 0.4 0.5 



FIG. 6. Results from fitting the model to the experimental data: non-racemic mixtures of Cis — TE — Glu. The solid line 
represents the linear correlation between experimental and numerical data obtained from fitting. The dotted line has slope 
equal to unity. 



TABLE VI. Results for the copolymerization model and experimental data correlations for non-racemic mixtures of C\% — 
TE - Lys. 



r 


Copolymerization model 


Binomial 


Fits for each subgroup n 


Global 
fit 


Global 
fit 


di 


tri 


tetra 


penta 


hexa 


hepta 


(R:S) 4:6 


0.78 


1 


0.87 


0.90 


0.84 


0.97 


0.89 




(R:S) 3:7 


0.93 


1 


0.95 


0.97 


0.99 




0.94 


0.65 



system polymerizes randomly. In Figure[9]the relative abundances of the oligopeptides are shown. The data correlation 
is shown in FigflOl 

The best fits for both chirally enriched mixture cases (20% and 40% excesses) in the case of Cx% — Glu — NCA are 
found for the same dynamics, that is k aa — — 3.8(s~ 1 mol~ 1 ). The results for these values are shown in Table 
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TABLE VII. Enantiomeric excesses: numerical results from the copolymerization model (experimental data) for the relative 
abundances of the homochiral oligopeptides for Cis — TE — Lys. 



ee(%) 


di 


tri 


tetra 


penta 


hexa 


hepta 


(R:S) 4:6 


23 (34) 


30 (34) 


36 (41) 


42 (60) 


49 (62) 


54 (>99.8) 


(R:S) 3:7 


45 (46) 


57 (63) 


66 (73) 


75 (85) 


81 (86) 






FIG. 7. Relative abundance versus number of repeat units (r, s) of the oligopeptides obtained from fitting the model (white) 
to the experimental data (black) from non-racemic mixtures of monomers of Ci8 — TE — Lys. 




FIG. 8. Results from fitting the model to the experimental data. Chiral mixtures of Cis — TE — Lys. The solid line represents 
the linear correlation between experimental and numerical data obtained from fitting. The dotted line has slope equal to unity. 



IIXI In FigQT]we compare the best fit against the experimentally obtained relative abundances of the oligopeptides. 
The corresponding data correlation is shown in Fig ll2l Numerical results for the racemic case have been found for 
< = 10 n (s). 
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TABLE VIII. Results for the copolymerization model and experimental data correlations for non-racemic mixtures of Cis — 
thio — Glu. 



r 


Copolymerization model 


Binomial 


Fits for each subgroup n 


Global 
fit 


Global 
fit 


di 


tri 


tetra 


penta 


hexa 


hepta 


(R:S) 4:6 


0.93 


0.98 


0.93 


0.92 


0.92 


0.91 


0.91 


0.93 


(R:S) 3:7 


0.89 


1 


0.99 


0.99 


0.98 




0.96 


0.97 




FIG. 9. Relative abundances versus number of repeat units (r, s) of the oligopeptides obtained from fitting the model (white) 
to the experimental data (black) for the non-racemic mixtures of Cis — thio ~ Glu. 



Data Correlation 




C18-thio-Glu4:6 
y = x 

y = 0.892 x 
r = 0.91 



C18-thio-Glu3:7 
y = x 
-y = 0.899 x 
= 0.95 



0.2 0.3 0.4 

Model 



FIG. 10. Results from fitting the model to the experimental data for non-racemic mixtures of the Cis — thio — Glu system. 
The solid line represents the linear correlation between experimental and numerical data obtained from fitting. The dotted line 
has slope equal to unity. 



IV. CONCLUSIONS 



The overall scheme for the chiral amplification process leading to the experimental data investigated here involves a 
self-assembly step followed by a lattice-controlled polymerization 0, [H| . It is this subsequent polymerization which 
is the prime focus of this paper. The authors of the experimental work stress that it is not at all straightforward 
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TABLE IX. Results for the copolymerization model and experimental data correlations for Cis — Glu — NCA. The global fit 
from the binomial distribution is shown for comparison. 



r 


Copolymerization model 


Binomial 


Fits for each subgroup n 


Global 
fit 


Global 
fit 


di 


tri 


tetra 


penta 


hexa 


hepta 


octa 


nona 


deca 


(R:S) 4:6 


-0.79 


0.9 


0.63 


0.74 


0.95 


0.89 


0.77 


0.86 




0.68 


0.11 


(R:S) 3:7 


-0.33 


0.79 


0.81 


0.76 


0.86 


0.96 


0.75 


0.83 


0.89 


0.75 






FIG. 11. Relative abundances for the non racemic mixtures of Ci&—thio — Glu. The experimental data set (black) is incomplete, 
we have used our model to fill in the missing portions of the histogram (white) 




Model 



FIG. 12. Results from fitting the model to the experimental data for the Cis — Glu — NCA system. The solid line represents 
the linear correlation between experimental and numerical data obtained from fitting. The dotted line has slope equal to unity. 
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to actually establish the correlation between the packing arrangement of the crystallites and the composition of the 
diastereoisomeric products that result therefrom. Therefore, our task here was to fit the outcome of these latter steps 
assuming an effective copolymerization scheme. The term "effective" simply means that the putative complicated 
correlations and interplay between the 2D crystallite phases at the air- water interface and the polymerization reaction 
pathways that depend on the microscopic packing arrangements within the crystals are treated here with a simple 
model. In this regard, our model can be regarded as a "course-grained" description of the overall process in that the 
microscopic details (the structures of the crystalline phases) are not resolved, but that the end-result or net effect of 
the pathways afforded by the crystallites can be summarized by the polymerization scheme as depicted graphically in 

Fig. m 

The model as introduced is defined for fully reversible reactions and this implies that some of the reaction rates must 
obey a corresponding constraint as dictated by microreversibility. Thus the model is appropriate for closed systems 
under thermodynamic control. For the numerical fits themselves, we found that all the reverse reaction rates could be 
set to rather tiny values, and this in consonance with experimentally observed irreversible condensation. Thus for the 
present purposes, the copolymerization model is practically irreversible. The values for the forward rates of adding 
the same chirality monomer to the end of the growing chain are found to be greater than those for addition of a wrong 
chirality monomer: that is, k aa = kbb > k a b = kt, a = 1 (except of course for the model system C\% — thio— Glu serving 
as reference for random polymerization). 

Other closed systems that lead to copolymers could be in principle be fit with our model. If for example k aa and 
kbb had different magnitudes, this would imply that a underlying chiral bias is operative either in the polymerization 
or in the prior formation of the two crystallites that control the polymerization. This bias could affect the packing 
arrangements of the crystal monomers and the reaction pathways taken within each crystallite phase. Since however 
our model is effective, as explained earlier, we would not be able to say whether the chiral bias is in the polymerization 
or in the structure of the crystallites that control the polymerization. Nevertheless, this bias in k aa being different 
from kbb, would result in favoring the attachment of say, an S to an S over the attachment of an R to an R, and this 
feature would show clearly up in the relative abundances. 

Another positive feature of the model is the robustness of the fits with respect to differing initial imbalances of 
the enantiomers. That is, for a given chemical model (including catalyst, if any) the values of the fitted rates do 
not depend on the initial enantiomeric excesses of the monomers. If our rate constants are viewed as effective, that 
is, implicitly involving the different chemical properties of the racemic and enantiomorphous crystallite phases, then 
this feature suggests that the packing arrangements and reaction pathways in the solid-state do not depend (or only 
weakly) on the magnitude of these imbalances. 

The Pearson product-moment correlation coefficient r between experimental and numerical data is greater for the 
copolymerization model than for the binomial distribution, except for the C\% — thio — Glu, which truly polymerizes 
randomly. The correlation between calculated and experimental relative abundances is also greater for the initially 
non-racemic situations, and the higher the initially chiral enrichment of the mixture is, the better the copolymerization 
model reproduces the chemical data. The results obtained here lead us to affirm that the model systems considered 
all undergo a non-random polymerization, as was asserted by the authors of the experiments [24], [25| • 

The model also qualitatively reproduces the behavior of the enantiomeric excess ee, its increase with the length of 
the chains and the enhancement of the ee of the corresponding initial mixture of monomers. All this, in spite of the 
complexity of the factors that affect the reactivity within the experimental two-phase system, i.e., the microscopic 
crystallite packing arrangements and the possible reaction pathways within these 2D crystallites. In conclusion then, 
we may therefore assert that our simple scheme does provide an accurate course-grained description of the lattice- 
controlled polymerization reported in Ref [24| HH ■ 



ACKNOWLEDGEMENTS 

We are grateful to Meir Lahav for providing us with the experimental data and for many helpful discussions and 
correspondence. CB has a Calvo-Rodes predoctoral scholarship from the Instituto Nacional de Tecnica Aeroespacial 
(INTA) and the research of DH is supported in part by the grant AYA2009-13920-C02-01 from the Ministerio de 
Ciencia e Innovation (Spain) and forms part of the COST Action CM0703 "Systems Chemistry" . 



[1] G. Joyce, G. Visser, C. van Boeckel, J. van Boom, L. Orgel, and J. van Westrenen, Nature 310, 602 (1984). 
[2] V. Avetisov and V. Goldanskii, Proceedings of the National Academy of Sciences of the United States of America 93, 
11435 (1996). 



18 



[3] D. Kondepudi and K. Asakura, Acc. Chem. Res. 34, 946 (2001). 
[4] P. Cintas, Angew. Chem. Int. Ed. 41, 1139 (2002). 

[5] C. Huber, W. Eisenreich, S. Hecht, and G. Wachtershauser, Science 301, 938 (2003). 
[6] L. Leman, L. Orgel, and M. Ghadiri, Science 306, 283 (2004). 
[7] B. 

[8] R. Pascal, L. Boiteau, and A. Commeyras, Top. Curr. Chem. 259, 69 (2005). 

[9] R. Lundberg and P. Doty, J. Am. Chem. Soc. 79, 3961 (1957). 
[10] M. Idelson and E. Blout, J. Am. Chem. Soc. 79, 3948 (1957). 
[11] T. Akaike and S. Inoue, Biopolymers 15, 1863 (1976). 
[12] N. Blair and W. Bonner, Orig. Life Evol. Biosph. 10, 255 (1980). 
[13] N. Blair and W. Bonner, Orig. Life Evol. Biosph. 11, 331 (1981). 
[14] H. Kricheldorf, Angew. Chem. Int. Ed. 45, 5752 (2006). 
[15] K. Ehler and L. Orgel, Biochem. et Biosphy. Acta 434, 233 (1976). 
[16] A. Hill, C. Bohler, and L. Orgel, Orig. Life Evol. Biosph. 28, 235 (1998). 
[17] A. Brack, Orig. Life Evol. Biosph. 17, 367 (1987). 
[18] H. Kanazawa, Polymer 33, 2557 (1992). 

[19] H. Kanazawa and Y. Ohashi, Mol. Cryst. Liq. Cryst. 277, 45 (1996). 

[20] T. Hitz, M. Blocher, P. Walde, and P. Luisi, Macromolecules 34, 2443 (2001). 

[21] T. Hitz and P. Luisi, Helv. Chim. Acta 85, 3975 (2002). 

[22] T. Hitz and P. Luisi, Helvetica Chimica Acta 86, 1423 (2003). 

[23] M. Blocher, T. Hitz, and P. Luisi, Helvetica Chimica Acta 84, 842 (2001). 

[24] H. Zepik, E. Shavit, M. Tang, T. Jensen, K. Kjaer, G. Bolbach, L. Leiserowitz, I. Weissbuch, and M. Lahav, Science 295, 
1266 (2002). 

[25] I. Weissbuch, H. Zepik, G. Bolbach, E. Shavit, M. Tang, T. Jensen, K. Kjaer, L. Leiserowitz, and M. Lahav, Chemistry: 

a European Journal 9, 1782 (2003). 
[26] J. A. D. Wattis and P. V. Co veney, Journal of Physical Chemistry B 111, 9546 (2007). 
[27] C. Blanco and D. Hochberg, |Phys Chem Chem Phys 13, 839 (20"TT) 



